function [w,itilde]=funcEM2(const,h,tmax,winit)
w = zeros(3,tmax);   % Preallocate w
w_= zeros(3,tmax);   % Preallocate w_
w(:,1) = winit; % Initial values
bool=false;
itilde=1;


for i = 2:h:tmax
    w_(:,i) = w(:,i-1)+h*func2U82(w(:,i-1),const,const.P_in(i)); % Predictor
    w(:,i) = w(:,i-1)+(h/2)*(func2U82(w(:,i-1),const,const.P_in(i))+func2U82(w_(:,i),const,const.P_in(i))); % Corrector
    
    if (w(1,i)<=.005 && bool==false)
        itilde=i;
        bool=true;
    else
        if (w(1,i)>=.005 && bool==true)
            bool=false;
            itilde=-1;
        end
    end
end